Load the necessary libraries
library(car) # for regression diagnostics
library(broom) # for tidy output
library(ggfortify) # for model diagnostics
library(sjPlot) # for outputs
library(knitr) # for kable
library(effects) # for partial effects plots
library(ggeffects) # for partial effects plots
library(emmeans) # for estimating marginal means
library(MASS) # for glm.nb
library(MuMIn) # for AICc
library(tidyverse) # for data wrangling
library(broom.mixed)
library(nlme) # for lme
library(lme4) # for lmer
library(lmerTest) # for Satterthwaite's p-values
library(glmmTMB) # for glmmTMB
library(DHARMa) # for residuals and diagnostics
library(performance) # for diagnostic plots
library(see) # for diagnostic plots
Starlings
Sampling design
Format of starling_full.RSV data files
| SITUATION | MONTH | MASS | BIRD |
|---|---|---|---|
| tree | Nov | 78 | tree1 |
| .. | .. | .. | .. |
| nest-box | Nov | 78 | nest-box1 |
| .. | .. | .. | .. |
| inside | Nov | 79 | inside1 |
| .. | .. | .. | .. |
| other | Nov | 77 | other1 |
| .. | .. | .. | .. |
| tree | Jan | 85 | tree1 |
| .. | .. | .. | .. |
| SITUATION | Categorical listing of roosting situations (tree, nest-box, inside or other) |
| MONTH | Categorical listing of the month of sampling. |
| MASS | Mass (g) of starlings. |
| BIRD | Categorical listing of individual bird repeatedly sampled. |
This is a split-plot (or repeated measures) design. The individual birds are the blocks, the Situation is the between block effect and the Month is the within block effect. Repeated measures analyses involve a within block effect that represents time (in this case Month). Since it is not possible to randomise the order of time, repeated measures designs have the potential for the residuals to be auto-correlated. That is, rather than being independent, residuals from observations that are closer in time, tend to be more similar (correlated) than the residuals associated with observations that are further apart in time.
That said, with only two time points, auto-correlation is not possible.
starling <- read_csv("../public/data/starling_full.csv", trim_ws = TRUE)
glimpse(starling)
## Rows: 80
## Columns: 5
## $ MONTH <chr> "Nov", "Nov", "Nov", "Nov", "Nov", "Nov", "Nov", "Nov", "No…
## $ SITUATION <chr> "tree", "tree", "tree", "tree", "tree", "tree", "tree", "tr…
## $ subjectnum <dbl> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 1, 2, 3, 4, 5, 6, 7, 8, 9, 1…
## $ BIRD <chr> "tree1", "tree2", "tree3", "tree4", "tree5", "tree6", "tree…
## $ MASS <dbl> 78, 88, 87, 88, 83, 82, 81, 80, 80, 89, 78, 78, 85, 81, 78,…
Lets prepare the data:
starling <- starling %>% mutate(
MONTH = factor(MONTH, levels = c("Nov", "Jan")),
SITUATION = factor(SITUATION),
BIRD = factor(BIRD)
)
Model formula: \[ y_i \sim{} \mathcal{N}(\mu_i, \sigma^2)\\ \mu_i =\boldsymbol{\beta} \bf{X_i} + \boldsymbol{\gamma} \bf{Z_i} \]
where \(\boldsymbol{\beta}\) and \(\boldsymbol{\gamma}\) are vectors of the fixed and random effects parameters respectively and \(\bf{X}\) is the model matrix representing the overall intercept and effects of roosting situation and month on starling mass. \(\bf{Z}\) represents a cell means model matrix for the random intercepts associated with individual birds.
ggplot(starling, aes(y = MASS, x = MONTH)) +
geom_boxplot() +
facet_grid(~SITUATION)
## Better still
ggplot(starling, aes(y = MASS, x = MONTH, group = BIRD)) +
geom_point() +
geom_line() +
facet_grid(~SITUATION)
Conclusions:
starling.lme1 <- lme(MASS ~ 1, random = ~ 1 | BIRD, data = starling, method = "REML")
starling.lme2 <- lme(MASS ~ 1, random = ~ MONTH | BIRD, data = starling, method = "REML")
AICc(starling.lme1, starling.lme2)
Conclusions:
Although we will explore this, we might as well just fit the multiplicative model and then explore the interaction from the outputs.
## Multiplicative and additive models
starling.lme1a <- update(starling.lme1, . ~ . + MONTH * SITUATION, method = "ML")
starling.lme1b <- update(starling.lme1, . ~ . + MONTH + SITUATION, method = "ML")
AICc(starling.lme1a, starling.lme1b)
## Finally, run as REML
starling.lme1a <- update(starling.lme1a, method = "ML")
Conclusions:
starling.lmer1 <- lmer(MASS ~ 1 + (1 | BIRD), data = starling, REML = TRUE)
starling.lmer2 <- lmer(MASS ~ 1 + (MONTH | BIRD), data = starling, REML = TRUE, control = lmerControl(check.nobs.vs.nRE = "ignore"))
AICc(starling.lmer1, starling.lmer2)
Conclusions:
## Multiplicative and additive models
starling.lmer1a <- update(starling.lmer1, . ~ . + MONTH * SITUATION, REML = FALSE)
starling.lmer1b <- update(starling.lmer1, . ~ . + MONTH + SITUATION, REML = FALSE)
AICc(starling.lmer1a, starling.lmer1b)
## Finally, run as REML
starling.lmer1a <- update(starling.lmer1a, REML = TRUE)
Conclusions:
starling.glmmTMB1 <- glmmTMB(MASS ~ 1 + (1 | BIRD), data = starling, REML = TRUE)
starling.glmmTMB2 <- glmmTMB(MASS ~ 1 + (MONTH | BIRD), data = starling, REML = TRUE)
starling.glmmTMB2 <- glmmTMB(MASS ~ 1 + (MONTH | BIRD),
data = starling, REML = TRUE,
control = glmmTMBControl(
optimizer = optim,
optArgs = list(method = "BFGS")
)
)
AICc(starling.glmmTMB1, starling.glmmTMB2)
Conclusions:
## Multiplicative and additive models
starling.glmmTMB1a <- update(starling.glmmTMB1, . ~ . + MONTH * SITUATION, REML = FALSE)
starling.glmmTMB1b <- update(starling.glmmTMB1, . ~ . + MONTH + SITUATION, REML = FALSE)
AICc(starling.glmmTMB1a, starling.glmmTMB1b)
## Finally, run as REML
starling.glmmTMB1a <- update(starling.glmmTMB1a, REML = TRUE)
Conclusions:
plot_model(starling.lme1a, type = "diag")[-2] %>% plot_grid()
Conclusions:
starling.lme1a %>% performance::check_model()
## Could not compute standard errors from random effects for diagnostic plot.
## Homogeneity of variance could not be computed. Cannot extract residual variance from objects of class 'lme'.
Conclusions:
# starling.resid = simulateResiduals(starling.lme1a, plot=TRUE)
plot_model(starling.lmer1a, type = "diag")[-2] %>% plot_grid()
Conclusions:
starling.lmer1a %>% performance::check_model()
Conclusions:
starling.resid <- starling.lmer1a %>% simulateResiduals(plot = TRUE)
Conclusions:
plot_model(starling.glmmTMB1a, type = "diag")[-2] %>% plot_grid()
Conclusions:
starling.glmmTMB1a %>% performance::check_model()
Conclusions:
starling.resid <- starling.glmmTMB1a %>% simulateResiduals(plot = TRUE)
Conclusions:
starling.lme1a %>% plot_model(type = "eff", terms = c("SITUATION", "MONTH"))
starling.lme1a %>% plot_model(type = "est")
starling.lme1a %>%
allEffects() %>%
plot(multiline = TRUE, ci.style = "bars")
starling.lme1a %>%
ggpredict(terms = c("SITUATION", "MONTH")) %>%
plot(add.data = TRUE)
starling.lme1a %>%
ggemmeans(~ SITUATION * MONTH) %>%
plot(add.data = TRUE)
starling.lmer1a %>% plot_model(type = "eff", terms = c("SITUATION", "MONTH"), show.data = TRUE)
starling.lmer1a %>% plot_model(type = "est")
starling.lmer1a %>%
allEffects() %>%
plot(multiline = TRUE, ci.style = "bars")
starling.lmer1a %>%
ggpredict(terms = c("SITUATION", "MONTH")) %>%
plot(add.data = TRUE)
starling.lmer1a %>%
ggemmeans(~ SITUATION * MONTH) %>%
plot(add.data = TRUE)
starling.glmmTMB1a %>% plot_model(type = "eff", terms = c("SITUATION", "MONTH"), show.data = TRUE)
starling.glmmTMB1a %>% plot_model(type = "est")
starling.glmmTMB1a %>%
allEffects() %>%
plot(multiline = TRUE, ci.style = "bars")
starling.glmmTMB1a %>%
ggpredict(terms = c("SITUATION", "MONTH")) %>%
plot(add.data = TRUE)
starling.glmmTMB1a %>%
ggemmeans(~ SITUATION * MONTH) %>%
plot(add.data = TRUE)
starling.lme1a %>% summary()
## Linear mixed-effects model fit by maximum likelihood
## Data: starling
## AIC BIC logLik
## 468.4462 492.2665 -224.2231
##
## Random effects:
## Formula: ~1 | BIRD
## (Intercept) Residual
## StdDev: 0.5567771 3.951582
##
## Fixed effects: MASS ~ MONTH + SITUATION + MONTH:SITUATION
## Value Std.Error DF t-value p-value
## (Intercept) 78.6 1.330205 36 59.08865 0.0000
## MONTHJan 9.6 1.862794 36 5.15355 0.0000
## SITUATIONnest-box 0.8 1.881193 36 0.42526 0.6732
## SITUATIONother -3.2 1.881193 36 -1.70105 0.0976
## SITUATIONtree 5.0 1.881193 36 2.65789 0.0117
## MONTHJan:SITUATIONnest-box 1.2 2.634388 36 0.45551 0.6515
## MONTHJan:SITUATIONother -0.8 2.634388 36 -0.30368 0.7631
## MONTHJan:SITUATIONtree -2.4 2.634388 36 -0.91103 0.3683
## Correlation:
## (Intr) MONTHJn SITUATION- SITUATIONth SITUATIONtr
## MONTHJan -0.700
## SITUATIONnest-box -0.707 0.495
## SITUATIONother -0.707 0.495 0.500
## SITUATIONtree -0.707 0.495 0.500 0.500
## MONTHJan:SITUATIONnest-box 0.495 -0.707 -0.700 -0.350 -0.350
## MONTHJan:SITUATIONother 0.495 -0.707 -0.350 -0.700 -0.350
## MONTHJan:SITUATIONtree 0.495 -0.707 -0.350 -0.350 -0.700
## MONTHJ:SITUATION- MONTHJn:SITUATIONth
## MONTHJan
## SITUATIONnest-box
## SITUATIONother
## SITUATIONtree
## MONTHJan:SITUATIONnest-box
## MONTHJan:SITUATIONother 0.500
## MONTHJan:SITUATIONtree 0.500 0.500
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -1.85043993 -0.81028549 -0.09107779 0.74016509 2.28662515
##
## Number of Observations: 80
## Number of Groups: 40
Conclusions:
starling.lme1a %>% fixef()
## (Intercept) MONTHJan
## 78.6 9.6
## SITUATIONnest-box SITUATIONother
## 0.8 -3.2
## SITUATIONtree MONTHJan:SITUATIONnest-box
## 5.0 1.2
## MONTHJan:SITUATIONother MONTHJan:SITUATIONtree
## -0.8 -2.4
starling.lme1a %>% tidy(conf.int = TRUE)
starling.lme1a %>%
tidy(conf.int = TRUE) %>%
kable()
| effect | group | term | estimate | std.error | df | statistic | p.value | conf.low | conf.high |
|---|---|---|---|---|---|---|---|---|---|
| fixed | fixed | (Intercept) | 78.6000000 | 1.330205 | 36 | 59.0886517 | 0.0000000 | 76.0406611 | 81.1593389 |
| fixed | fixed | MONTHJan | 9.6000000 | 1.862794 | 36 | 5.1535501 | 0.0000094 | 6.0159500 | 13.1840500 |
| fixed | fixed | SITUATIONnest-box | 0.8000000 | 1.881193 | 36 | 0.4252619 | 0.6731771 | -2.8194518 | 4.4194518 |
| fixed | fixed | SITUATIONother | -3.2000000 | 1.881193 | 36 | -1.7010478 | 0.0975578 | -6.8194518 | 0.4194518 |
| fixed | fixed | SITUATIONtree | 5.0000000 | 1.881193 | 36 | 2.6578872 | 0.0116541 | 1.3805482 | 8.6194518 |
| fixed | fixed | MONTHJan:SITUATIONnest-box | 1.2000000 | 2.634388 | 36 | 0.4555138 | 0.6514753 | -3.8686122 | 6.2686122 |
| fixed | fixed | MONTHJan:SITUATIONother | -0.8000000 | 2.634388 | 36 | -0.3036759 | 0.7631231 | -5.8686122 | 4.2686122 |
| fixed | fixed | MONTHJan:SITUATIONtree | -2.4000000 | 2.634388 | 36 | -0.9110276 | 0.3683410 | -7.4686122 | 2.6686122 |
| ran_pars | BIRD | sd_(Intercept) | 0.5567771 | NA | NA | NA | NA | 0.0000881 | 3520.1759737 |
| ran_pars | Residual | sd_Observation | 3.9515819 | NA | NA | NA | NA | NA | NA |
Conclusions:
# warning this is only appropriate for html output
starling.lme1a %>% sjPlot::tab_model(show.se = TRUE, show.aic = TRUE)
| Â | MASS | |||
|---|---|---|---|---|
| Predictors | Estimates | std. Error | CI | p |
| (Intercept) | 78.60 | 1.33 | 76.04 – 81.16 | <0.001 |
| MONTH [Jan] | 9.60 | 1.86 | 6.02 – 13.18 | <0.001 |
| SITUATION [nest-box] | 0.80 | 1.88 | -2.82 – 4.42 | 0.673 |
| SITUATION [other] | -3.20 | 1.88 | -6.82 – 0.42 | 0.098 |
| SITUATION [tree] | 5.00 | 1.88 | 1.38 – 8.62 | 0.012 |
|
MONTH [Jan] * SITUATION [nest-box] |
1.20 | 2.63 | -3.87 – 6.27 | 0.651 |
|
MONTH [Jan] * SITUATION [other] |
-0.80 | 2.63 | -5.87 – 4.27 | 0.763 |
|
MONTH [Jan] * SITUATION [tree] |
-2.40 | 2.63 | -7.47 – 2.67 | 0.368 |
| Random Effects | ||||
| σ2 | 15.61 | |||
| τ00 BIRD | 0.31 | |||
| ICC | 0.02 | |||
| N BIRD | 40 | |||
| Observations | 80 | |||
| Marginal R2 / Conditional R2 | 0.643 / 0.650 | |||
| AIC | 468.446 | |||
starling.lmer1a %>% summary()
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: MASS ~ (1 | BIRD) + MONTH + SITUATION + MONTH:SITUATION
## Data: starling
##
## REML criterion at convergence: 429.6
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.7555 -0.7687 -0.0864 0.7022 2.1693
##
## Random effects:
## Groups Name Variance Std.Dev.
## BIRD (Intercept) 0.3444 0.5869
## Residual 17.3500 4.1653
## Number of obs: 80, groups: BIRD, 40
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 78.600 1.330 71.973 59.089 < 2e-16 ***
## MONTHJan 9.600 1.863 36.000 5.154 9.39e-06 ***
## SITUATIONnest-box 0.800 1.881 71.973 0.425 0.67191
## SITUATIONother -3.200 1.881 71.973 -1.701 0.09325 .
## SITUATIONtree 5.000 1.881 71.973 2.658 0.00968 **
## MONTHJan:SITUATIONnest-box 1.200 2.634 36.000 0.456 0.65148
## MONTHJan:SITUATIONother -0.800 2.634 36.000 -0.304 0.76312
## MONTHJan:SITUATIONtree -2.400 2.634 36.000 -0.911 0.36834
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) MONTHJn SITUATION- SITUATIONth SITUATIONtr
## MONTHJan -0.700
## SITUATIONn- -0.707 0.495
## SITUATIONth -0.707 0.495 0.500
## SITUATIONtr -0.707 0.495 0.500 0.500
## MONTHJ:SITUATION- 0.495 -0.707 -0.700 -0.350 -0.350
## MONTHJn:SITUATIONth 0.495 -0.707 -0.350 -0.700 -0.350
## MONTHJn:SITUATIONtr 0.495 -0.707 -0.350 -0.350 -0.700
## MONTHJ:SITUATION- MONTHJn:SITUATIONth
## MONTHJan
## SITUATIONn-
## SITUATIONth
## SITUATIONtr
## MONTHJ:SITUATION-
## MONTHJn:SITUATIONth 0.500
## MONTHJn:SITUATIONtr 0.500 0.500
Conclusions:
confint(starling.lmer1a)
## 2.5 % 97.5 %
## .sig01 0.000000 2.4290757
## .sigma 3.221630 4.6973017
## (Intercept) 76.096636 81.1033641
## MONTHJan 6.066732 13.1332685
## SITUATIONnest-box -2.740292 4.3402915
## SITUATIONother -6.740292 0.3402915
## SITUATIONtree 1.459708 8.5402915
## MONTHJan:SITUATIONnest-box -3.796795 6.1967962
## MONTHJan:SITUATIONother -5.796795 4.1967962
## MONTHJan:SITUATIONtree -7.396795 2.5967962
starling.lmer1a %>% tidy(conf.int = TRUE)
starling.lmer1a %>%
tidy(conf.int = TRUE) %>%
kable()
| effect | group | term | estimate | std.error | statistic | df | p.value | conf.low | conf.high |
|---|---|---|---|---|---|---|---|---|---|
| fixed | NA | (Intercept) | 78.6000000 | 1.330205 | 59.0886517 | 71.97273 | 0.0000000 | 75.992847 | 81.2071532 |
| fixed | NA | MONTHJan | 9.6000000 | 1.862794 | 5.1535500 | 35.99999 | 0.0000094 | 5.948992 | 13.2510084 |
| fixed | NA | SITUATIONnest-box | 0.8000000 | 1.881193 | 0.4252619 | 71.97273 | 0.6719146 | -2.887071 | 4.4870715 |
| fixed | NA | SITUATIONother | -3.2000000 | 1.881193 | -1.7010478 | 71.97273 | 0.0932510 | -6.887072 | 0.4870715 |
| fixed | NA | SITUATIONtree | 5.0000000 | 1.881193 | 2.6578872 | 71.97273 | 0.0096818 | 1.312928 | 8.6870715 |
| fixed | NA | MONTHJan:SITUATIONnest-box | 1.2000000 | 2.634388 | 0.4555138 | 35.99999 | 0.6514753 | -3.963306 | 6.3633056 |
| fixed | NA | MONTHJan:SITUATIONother | -0.8000000 | 2.634388 | -0.3036758 | 35.99999 | 0.7631231 | -5.963306 | 4.3633056 |
| fixed | NA | MONTHJan:SITUATIONtree | -2.4000000 | 2.634388 | -0.9110275 | 35.99999 | 0.3683410 | -7.563306 | 2.7633056 |
| ran_pars | BIRD | sd__(Intercept) | 0.5868937 | NA | NA | NA | NA | NA | NA |
| ran_pars | Residual | sd__Observation | 4.1653331 | NA | NA | NA | NA | NA | NA |
Conclusions:
# warning this is only appropriate for html output
starling.lmer1a %>% sjPlot::tab_model(show.se = TRUE, show.aic = TRUE)
| Â | MASS | |||
|---|---|---|---|---|
| Predictors | Estimates | std. Error | CI | p |
| (Intercept) | 78.60 | 1.33 | 75.99 – 81.21 | <0.001 |
| MONTH [Jan] | 9.60 | 1.86 | 5.95 – 13.25 | <0.001 |
| SITUATION [nest-box] | 0.80 | 1.88 | -2.89 – 4.49 | 0.671 |
| SITUATION [other] | -3.20 | 1.88 | -6.89 – 0.49 | 0.089 |
| SITUATION [tree] | 5.00 | 1.88 | 1.31 – 8.69 | 0.008 |
|
MONTH [Jan] * SITUATION [nest-box] |
1.20 | 2.63 | -3.96 – 6.36 | 0.649 |
|
MONTH [Jan] * SITUATION [other] |
-0.80 | 2.63 | -5.96 – 4.36 | 0.761 |
|
MONTH [Jan] * SITUATION [tree] |
-2.40 | 2.63 | -7.56 – 2.76 | 0.362 |
| Random Effects | ||||
| σ2 | 17.35 | |||
| τ00 BIRD | 0.34 | |||
| ICC | 0.02 | |||
| N BIRD | 40 | |||
| Observations | 80 | |||
| Marginal R2 / Conditional R2 | 0.618 / 0.626 | |||
| AIC | 449.608 | |||
starling.glmmTMB1a %>% summary()
## Family: gaussian ( identity )
## Formula: MASS ~ (1 | BIRD) + MONTH + SITUATION + MONTH:SITUATION
## Data: starling
##
## AIC BIC logLik deviance df.resid
## 449.6 473.4 -214.8 429.6 78
##
## Random effects:
##
## Conditional model:
## Groups Name Variance Std.Dev.
## BIRD (Intercept) 0.3444 0.5869
## Residual 17.3500 4.1653
## Number of obs: 80, groups: BIRD, 40
##
## Dispersion estimate for gaussian family (sigma^2): 17.3
##
## Conditional model:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 78.600 1.330 59.09 < 2e-16 ***
## MONTHJan 9.600 1.863 5.15 2.56e-07 ***
## SITUATIONnest-box 0.800 1.881 0.43 0.67065
## SITUATIONother -3.200 1.881 -1.70 0.08893 .
## SITUATIONtree 5.000 1.881 2.66 0.00786 **
## MONTHJan:SITUATIONnest-box 1.200 2.634 0.46 0.64874
## MONTHJan:SITUATIONother -0.800 2.634 -0.30 0.76137
## MONTHJan:SITUATIONtree -2.400 2.634 -0.91 0.36228
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Conclusions:
starling.glmmTMB1a %>% confint()
## 2.5 % 97.5 % Estimate
## cond.(Intercept) 75.9928467763 81.2071532 78.6000000
## cond.MONTHJan 5.9489916440 13.2510084 9.6000000
## cond.SITUATIONnest-box -2.8870714481 4.4870714 0.8000000
## cond.SITUATIONother -6.8870714481 0.4870714 -3.2000000
## cond.SITUATIONtree 1.3129285519 8.6870714 5.0000000
## cond.MONTHJan:SITUATIONnest-box -3.9633055333 6.3633055 1.2000000
## cond.MONTHJan:SITUATIONother -5.9633055333 4.3633055 -0.8000000
## cond.MONTHJan:SITUATIONtree -7.5633055333 2.7633055 -2.4000000
## BIRD.cond.Std.Dev.(Intercept) 0.0001330334 2589.1574686 0.5868939
starling.glmmTMB1a %>% tidy(conf.int = TRUE)
starling.glmmTMB1a %>%
tidy(conf.int = TRUE) %>%
kable()
| effect | component | group | term | estimate | std.error | statistic | p.value | conf.low | conf.high |
|---|---|---|---|---|---|---|---|---|---|
| fixed | cond | NA | (Intercept) | 78.6000000 | 1.330205 | 59.0886519 | 0.0000000 | 75.992847 | 81.2071532 |
| fixed | cond | NA | MONTHJan | 9.6000000 | 1.862794 | 5.1535500 | 0.0000003 | 5.948992 | 13.2510084 |
| fixed | cond | NA | SITUATIONnest-box | 0.8000000 | 1.881193 | 0.4252619 | 0.6706457 | -2.887071 | 4.4870714 |
| fixed | cond | NA | SITUATIONother | -3.2000000 | 1.881193 | -1.7010478 | 0.0889340 | -6.887071 | 0.4870714 |
| fixed | cond | NA | SITUATIONtree | 5.0000000 | 1.881193 | 2.6578872 | 0.0078632 | 1.312929 | 8.6870714 |
| fixed | cond | NA | MONTHJan:SITUATIONnest-box | 1.2000000 | 2.634388 | 0.4555138 | 0.6487397 | -3.963306 | 6.3633055 |
| fixed | cond | NA | MONTHJan:SITUATIONother | -0.8000000 | 2.634388 | -0.3036758 | 0.7613749 | -5.963305 | 4.3633055 |
| fixed | cond | NA | MONTHJan:SITUATIONtree | -2.4000000 | 2.634388 | -0.9110275 | 0.3622809 | -7.563306 | 2.7633055 |
| ran_pars | cond | BIRD | sd__(Intercept) | 0.5868939 | NA | NA | NA | 0.000133 | 2589.1574686 |
| ran_pars | cond | Residual | sd__Observation | 4.1653331 | NA | NA | NA | 0.000133 | 2589.1574686 |
Conclusions:
# warning this is only appropriate for html output
starling.glmmTMB1a %>% sjPlot::tab_model(show.se = TRUE, show.aic = TRUE)
| Â | MASS | |||
|---|---|---|---|---|
| Predictors | Estimates | std. Error | CI | p |
| (Intercept) | 78.60 | 1.33 | 75.99 – 81.21 | <0.001 |
| MONTH [Jan] | 9.60 | 1.86 | 5.95 – 13.25 | <0.001 |
| SITUATION [nest-box] | 0.80 | 1.88 | -2.89 – 4.49 | 0.671 |
| SITUATION [other] | -3.20 | 1.88 | -6.89 – 0.49 | 0.089 |
| SITUATION [tree] | 5.00 | 1.88 | 1.31 – 8.69 | 0.008 |
|
MONTH [Jan] * SITUATION [nest-box] |
1.20 | 2.63 | -3.96 – 6.36 | 0.649 |
|
MONTH [Jan] * SITUATION [other] |
-0.80 | 2.63 | -5.96 – 4.36 | 0.761 |
|
MONTH [Jan] * SITUATION [tree] |
-2.40 | 2.63 | -7.56 – 2.76 | 0.362 |
| Random Effects | ||||
| σ2 | 17.35 | |||
| τ00 BIRD | 0.34 | |||
| ICC | 0.02 | |||
| N BIRD | 40 | |||
| Observations | 80 | |||
| Marginal R2 / Conditional R2 | 0.618 / 0.626 | |||
| AIC | 449.608 | |||
starling.lme1a %>%
emmeans(~SITUATION) %>%
pairs() %>%
summary(infer = TRUE)
starling.lme1a %>% r.squaredGLMM()
## R2m R2c
## [1,] 0.642884 0.6498357
## # R2 for Mixed Models
##
## Conditional R2: 0.650
## Marginal R2: 0.643
starling.lmer1a %>%
emmeans(~SITUATION) %>%
pairs() %>%
summary(infer = TRUE)
starling.lmer1a %>% r.squaredGLMM()
## R2m R2c
## [1,] 0.6183482 0.6257775
## # R2 for Mixed Models
##
## Conditional R2: 0.626
## Marginal R2: 0.618
starling.glmmTMB1a %>%
emmeans(~SITUATION) %>%
pairs() %>%
summary(infer = TRUE)
starling.glmmTMB1a %>% r.squaredGLMM()
## R2m R2c
## [1,] 0.6183482 0.6257776
## # R2 for Mixed Models
##
## Conditional R2: 0.626
## Marginal R2: 0.618
newdata <- starling.lme1a %>%
emmeans(~ SITUATION * MONTH) %>%
as.data.frame()
ggplot(newdata, aes(y = emmean, x = SITUATION, fill = MONTH)) +
geom_pointrange(aes(ymin = lower.CL, ymax = upper.CL),
shape = 21,
position = position_dodge(0.2)
) +
theme_classic()
newdata <- starling.lmer1a %>%
emmeans(~ SITUATION * MONTH) %>%
as.data.frame()
ggplot(newdata, aes(y = emmean, x = SITUATION, fill = MONTH)) +
geom_pointrange(aes(ymin = lower.CL, ymax = upper.CL),
shape = 21,
position = position_dodge(0.2)
) +
theme_classic()
newdata <- starling.glmmTMB1a %>%
emmeans(~ SITUATION * MONTH) %>%
as.data.frame()
ggplot(newdata, aes(y = emmean, x = SITUATION, fill = MONTH)) +
geom_pointrange(aes(ymin = lower.CL, ymax = upper.CL),
shape = 21,
position = position_dodge(0.2)
) +
theme_classic()